import sympy as sy
import numpy as np
import matplotlib.pyplot as plt

# EULER:
#for i in range(N):
#   aceler = g    -> queda livre
#   x[i+1] = x[i] + vx[i] * dt
#   vx[i+1] = vx[i] + aceler * dt
#   t[i+1] = t[i] + dt


#Considere a queda livre, sem resistência do ar.
#posição do objeto, usando o método de Euler [y = x, vy = vx pois e uma queda]

g = 9.8
dt = 0.1 #passo

t0 = 0
tf = 4

n = int((tf-t0)/dt) + 1 #soma 1 para incluir t = tf

t = np.linspace(0, 4, n)
y = np.empty(n)
vy = np.empty(n)

vy[0] = 0
t[0] = 0

for i in range(n-1):
    y[i + 1] = y[i] + vy[i] * dt
    vy[i + 1] = vy[i] + g * dt
    t[i + 1] = t[i] + dt


for i in range(n):
    print(f"t = {t[i]:.2f} s, y = {y[i]:.2f} m, vy = {vy[i]:.2f} m/s")

print(f"Velocidade no instante t = 3s (dt=0.1): {vy[30]:.2f} m/s")


#Repita as alíneas anteriores, com um passo de tempo 10 vezes menor
print("\n")

dt = 0.01
n = int((tf-t0)/dt) + 1  #soma 1 para incluir t = tf

t = np.linspace(0, 4, n)
y = np.empty(n)
vy = np.empty(n)

vy[0] = 0

for i in range(n-1):

    y[i + 1] = y[i] + vy[i] * dt
    vy[i + 1] = vy[i] + g * dt
    t[i + 1] = t[i] + dt


for i in range(n):
    print(f"t = {t[i]:.2f} s, y = {y[i]:.2f} m, vy = {vy[i]:.2f} m/s")

print("\n")
print(f"Velocidade no instante t = 3s (dt=0.01): {vy[300]:.2f} m/s")



#Qual a posição em 3s, se o objeto partiu da posição 0 m?

y0 = 0
t0 = 0
tf = 4
dt = 0.1

n = int((tf-t0)/dt) + 1  #soma 1 para incluir t = tf

t = np.linspace(0, 4, n)
y = np.empty(n)
vy = np.empty(n)

y[0] = 0
vy[0] = 0

for i in range(n-1):

    y[i + 1] = y[i] + vy[i] * dt
    vy[i + 1] = vy[i] + g * dt
    t[i + 1] = t[i] + dt

print("\n")
print(f"Posicao no instante t = 3s, partindo de y = 0m (dt=0.1): {y[30]:.2f} m/s")



#Calcule novamente a posição no instante 2s.
#Faça o gráfico do desvio do valor aproximado com o valor exato em função do passo. 
#Como varia o erro com o passo?

t_target = 2.0 
y_exato = 0.5 * g * (t_target**2)  # posicao exata em t = 2s: 19.6 m

dt_values = [0.1, 0.01, 0.001]
erros = []


for dt in dt_values:
    t0 = 0
    tf = 4
    y0 = 0
    vy0 = 0

    n = int((tf-t0)/dt) + 1  #soma 1 para incluir t = tf

    t = np.linspace(0, 4, n)
    y = np.empty(n)
    vy = np.empty(n)

    y[0] = 0
    vy[0] = 0


    for i in range(n-1):
        
        y[i + 1] = y[i] + vy[i] * dt
        vy[i + 1] = vy[i] + g * dt
        t[i + 1] = t[i] + dt
    
    y_aprox = y[int(t_target / dt)]
    erro = abs(y_exato - y_aprox)
    erros.append(erro)
    print(f"dt = {dt:1.3f}s, y_exato = {y_exato:.1f}m, y_aprox = {y_aprox:.1f}m, erro = {erro:1f}m")



plt.figure()
plt.loglog(dt_values, erros, marker='o', linestyle='-')
plt.xlabel("Passo (dt) [s]")
plt.ylabel("Erro absoluto em y(2) [m]")
plt.title("Erro do metodo de Euler em funcao do passo de tempo")
plt.grid(True, which="both", ls="--")
plt.show()